# it is implemented with null model stuff
# 
# Author: Guochun Shen
# Data:   2011-12-30
# Project:spatial statistic
# Email:  shenguochun@gmail.com
###############################################################################

phylocorr_null=function(com,phyd,focal_mark=NULL,rmax,step=1,nsim=19,alpha=0.05,scale=TRUE){
	if(is.null(focal_mark)){
		focal=1:(com$n)-1
		nfocal=com$n
	}else{
		focal=which(get_marks(com,focal_mark))-1
		nfocal=length(focal)
	}
	sp=get_marks(com,marks="sp")
	tips=rownames(phyd)
	sp=match(sp,tips)-1
	tips=1:length(tips)-1
	ntotal=com$n
	nsp=length(tips)
	
	nsteps = ceiling((rmax)/(step));
	pk=rep(0,nsteps)
	upper=rep(0,nsteps)
	lower=rep(0,nsteps)
	mean=rep(0,nsteps)
	pvalues=rep(0,nsteps)
	
	scale=as.numeric(scale)
	
	re=.C("phylocorr2", as.integer(focal), as.integer(nfocal), as.double(com$x), as.double(com$y),
			as.integer(ntotal), as.integer(sp), as.integer(nsp), as.double(rmax), 
			as.double(as.numeric(phyd)), as.double(step), as.double(pk), as.integer(nsim),
			as.double(mean), as.double(upper), as.double(lower), as.double(pvalues),
			as.double(alpha), as.integer(scale));
	r=seq(step/2,rmax,step)
	result=list("real"=re[[11]][-1],"upper"=re[[14]][-1],"lower"=re[[15]][-1],"r"=r[-1],"pvalues"=re[[16]][-1])
	return(result)
}

